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^ . Abstract: samurai is a tool for the automated numerical evaluation of one- loop cor- 



rections to any scattering amplitudes within the dimensional-regularization scheme. It is 
based on the decomposition of the integrand according to the OPP-approach, extended to 
accommodate an implementation of the generalized d-dimensional unitarity-cuts technique, 
and uses a polynomial interpolation exploiting the Discrete Fourier Transform, samurai 
can process integrands written either as numerator of Feynman diagrams or as product 
of tree- level amplitudes. We discuss some applications, among which the 6- and 8-photon 
scattering in QED, and the 6-quark scattering in QCD. samurai has been implemented 
as a Fortran90 library, publicly available, and it could be a useful module for the sys- 
tematic evaluation of the virtual corrections oriented towards automating next-to-leading 
order calculations relevant for the LHC phenomenology. 
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1. Introduction 

With the beginning of the experimental programs at the LHC, the need of describing 
particle scattering events with high accuracy becomes more pressing. On the theoretical 
side, perturbative calculation within leading order precision cannot be sufficient, therefore 
accounting for effects due to next-to-leading order corrections becomes mandatory. 

Leading order (LO) processes are the core of well-established multi-purpose event 
generators hke MadGraph MadEvent [1-3], CompHEP-CalcHEP [4,5], SHERPA [6,7], 
WHIZARD [8], ALPGEN [9], and HELAC [10,11], whereas a variety of processes computed 
at NLO are currently implemented in programs like MCFM [12,13] and NLOJET++ [14], 
or MC@NLO [15, 16] and the POWHEG [17-23] which consider also the matching with 
parton showers. 

The next-to-leading order (NLO) corrections to an n-parton final state process receive 
contributions from two sources: the one-loop correction to the (2 — )• n)-scattering, due to 
the exchange of an internal virtual particle; and the tree-level scattering (2—7-71 + 1), due 
to the real emission of an extra parton. Each contribution contains divergencies which 
cancel mutually in the final result where they are combined. 

The extraction of the real radiation singularities from general processes has been ad- 
dressed with techniques based either on phase-space slicing [24] or on the use of universal 
subtraction terms [25,26], which at present have been implemented in several algorithms, 
like the FKS subtraction [27], dipole subtraction [28-36] and antenna subtraction [37-40]. 

The increasing computational complexity of one-loop amplitudes, when the number of 
particles involved in the scattering increases, has limited the possibility of developing an 
automated multi-process evaluator for scattering amplitudes at NLO. The available results 
have been so far computed on a process-by-process basis, but, due to the recent advances 
in computational techniques for high-energy physics, that possibility is now at the horizon. 

Currently, the state-of-the-art is represented by the numerical calculation of extremely 
challenging 2 — )• 4 processes, like the EW corrections to e^e~ — )• 4f [41,42], or the NLO 
QCD corrections to pp ^ W + 3 jet production [43-46], — )• Z + 3 jet production [47], 
pp — )• ttbb [48-50] and pp — )• tijj [51], and qq — )• bbbb [52], which have been obtained 
both by optimizing the algebraic tensor reduction, and by developing novel approaches 
based on properties of scattering amplitudes such as factorization and unitarity. Also, 
the development of novel analytic techniques has benefited from a more systematic use of 
unitarity-based methods in combination with the theory of multivariate-complex functions 
[53-62], yielding the recent completion of the one-loop QCD correction to — t- H + 2 jets 
in the heavy-top limit [63-69]. 

It is well known that any one-loop amplitude can be expressed as a linear combination 
of a limited set of Master Integrals (MI) [70, 71]: therefore, the evaluation of one- loop 
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corrections reduces to evaluating the coefficients that multiply each MI. Aiming at the full 
reconstruction of one-loop amplitudes through such a decomposition, several automated 
packages have appeared, either in public releases like CutTools [72] and Golem [73], or 
in private versions such as the routines described in [41], [74] and [75], and codes like 
BlackHat [76], Rocket [77], and Helac-lLoop [78]. 

The development of novel numerical techniques have received a boost by the combina- 
tion of three important ideas: 

i) universal four-dimensional decomposition for the numerator of the integrand for any 
one-loop scattering amplitudes [79,80]; 

a) four-dimensional unitarity-cuts, detecting only the (poly)logarithmic structure of the 
amplitude, known as the cut-constructible term [53, 81] (see [82] for a more compre- 
hensive list of references); 

iii) unitarity-cuts in d-dimension, yielding the complete determination of dimensionally 
regulated one- loop amplitudes [83-88]. 

The first two ideas merged in the by-now known as OPP-approach [80,89], proposed by 
Papadopoulos, Pittau, and one of us, where the multi-pole decomposition of the numerator 
of any Feynman integral is achieved by a polynomial sampling that exploits the solutions 
of generalized unitarity-cuts. 

In the context of four-dimensional unitarity, the problem of computing the cut con- 
structible term and the rational term, that escapes the four-dimensional detection, are 
necessarily considered as separate issues. The reconstruction of the latter usually requires 
information from an extra source. When not obtained from the direct calculation of Feyn- 
man integrals, the rational term can be reconstructed by adding a piece derived from the 
cut-constructibic part (for instance, the ovcrlapping-tcrm within the on-shcU method [90], 
or the i?i-tcrm within the OPP-approach [91]), and a remaining piece computed through 
an additional tree-level like construction (for instance, the BCFW-recursive term within 
the on-shell method [90], or the i?2-term within the OPP-approach [91-93]). 

The idea of performing unitarity-cuts in d-dimension [60, 83-88] yields a combined 
determination of both cut-constructible and rational terms at once. This technique has 
been neatly systematized for numerical purposes by Ellis, Giele, Kunszt and Melnikov 
[87,88], and later proposed also for the on-shell approach by Badger [60]. 

In this paper we present SAMURAI, a tool based on a hybrid algorithm for the nu- 
merical computation of one-loop amplitudes, samurai relies on the extension of the 
OPP-polynomial structures to include an explicit dependence on the extra-dimensional 
parameter needed for the automated computation of the full rational term according to 
the d-dimensional approach, and makes use of a polynomial interpolation based on the 
type of Discrete Fourier Transform (DFT) described in [94]. 

We aim at producing a versatile code which could deal with any one- loop corrections, 
in massless as well as massive theories. Our reduction algorithm can process both (numer- 
ator of) Feynman integrals, proper of diagrammatic methods, and products of tree-level 
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amplitudes, as adopted in the framework of unitarity-based techniques. 
For a complete reconstruction of the rational term, the input should contain an explicit 
dependence on the dimensional-regularization parameters. In fact, it is expected to have a 
polynomial behavior in fi^, being fx the radial integration variable in the extra-dimensional 
subspace, and in e ( = (4 — (i)/2 ) according to the choice of the regularization scheme. 
The result is given as Laurent expansion in e up to the finite-order, and accounts for the 
full rational terms. 

SAMURAI is implemented as a Fortran90 library, publicly available at the webpage: 



and it is linked to OneLOop [78] and QCDLoop [95] for the numerical evaluation of the MI. 
We applied it to a series of known processes, like the four-, six-photon and eight-photon 
scattering in QED, the QCD virtual corrections to Drell-Yan, to the leading-color am- 
plitude for V + Ijet production, to the six-quark scattering, qiqi — > g2?2 Q3Q3, and to the 
contributions of the massive-scalar loop-diagrams to the all-plus helicity five- and six-gluon 
scattering. 

In particular, for the virtual corrections to qiqi — t- 52^2 Q3Q3 [52], we also considered the 
reduction of automatically generated integrands, by interfacing samurai with an infras- 
tructure derived from golein-2.0 [96], which provides numerators of Feynman integrals. 
These examples are thought to be used both as a guide to understand the SAMURAI frame- 
work, and as templates to generate the codes for other calculations. 

In the context of collaborations among different groups aiming at automated NLO cal- 
culations relevant for LHC phenomenology [97], and, therefore, providing complementary 
structures to be interfaced [98], samurai could constitute the module for the systematic 
evaluation of the virtual corrections. 

The paper is organized as follows. The reduction algorithm is discussed in Section 2; 
Section 3 describes the key-points of the SAMURAI library, while a series of applications are 
illustrated in Section 4. In Section 5, we resume our conclusions. 

2. Reduction Algorithm 

The reduction method is based on the general decomposition for the integrand of a generic 
one-loop amplitude, originally proposed by Papadopoulos, Pittau and one of us [80,89], 
and later extended by Ellis, Giele, Kunszt and Melnikov [87,88]. Within the dimensional 
regularization scheme, any one-loop n-point amplitude can be written as 



We use a bar to denote objects living in d = 4 — 2e dimensions, following the prescription 
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Di = {q+ Pif - ml = {q + pif - mf - i?, 



{PO ^ 0) . 



(2.1) 



i= i + ii , with 




2 2 
q -11 . 



(2.2) 
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Also, we use the notation /(g) as short-hand notation for f{q,fJ-'^). 



2.1 Integrands 

SAMURAI can reduce integrands of one-loop amplitudes which can be defined in two ways, 
either as numerator functions (sitting on products of denominators), or as products of tree- 
level amplitudes (sewn along cut-lines). The former definition accommodates a reduction 
based on a diagrammatic method, while the latter is proper of a unitarity-based technology. 

According to the chosen dimensional regularization scheme, the most general numer- 
ator of one-loop amplitudes A/'(g, e) can be thought as composed of three terms, 

Af{q, e) = No{q) + eNi{q) + e^N2{q). (2.3) 

The coefficients of this e-expansion, Nq, Ni and are functions of q'^ and /i^, therefore 
in our discussion, except when a distinction between them is necessarily required, we will 
simply talk about N, giving as understood that the same logic would apply to each of the 
three contributions Ni. 

2.1.1 Decomposition 

According to [80,89], the numerator N{q) can be expressed in terms of denominators Di, 
as follows 

n— 1 n— 1 n— 1 n— 1 

N{q) = ^ijklUfi) n ^h+Y. ^^JkM n + 

i«m hyti,j,k,£,m i«£ hyti,j,k,£ 

n—l n—1 n—1 n— 1 n— 1 n— 1 

+ ^ A,^k{q) U Dh + Y, ^iM Wd^ + Y. ^^(^) n , (2.4) 

i«k hy^i,j,k i<j h^i,j i h^i 

where i « m stands for a lexicographic ordering i<j<k<i<m. The functions A(g) = 
A(g, /i^) are polynomials in the components of q and in fi^. By using the decomposition 
( |2.4| ) in Eq.( p.ip , the multi-pole nature of the integrand of any one-loop n-point amplitude 
becomes trivially exposed, 

(q) , sr^ \jki{q) , ^ijkiq) , 



A,,iq) , r.'A,;(g) 



i<j 'J i ' 

which, as we will see, is responsible of the decomposition of any dimensional regulated 
one- loop amplitude in terms of Master Integrals (MI) associated to 4-, 3-, 2-, and 1-point 
functions, respectively called boxes, triangles, bubbles, and tadpoles. 

2.2 Polynomial Structures and Discrete Fourier Transform 

The calculation of a generic scattering amplitude amounts to the problem of extracting 
the coefficients of multivariate polynomials, generated at every step of the multiple-cut 
analysis. To determine these coefficients we implement a semi-numerical algorithm whose 
main features are: 
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• the extension of the OPP-polynomials [80,89] for quadruple-, triple- and double-cut 
to the framework of d-dimensional unitarity [87, 88] ; 

• the parametrization of the residue of the quintuple-cut affecting only the polynomial 
dependence on the extra-dimension scale [99]; 

• the numerical sampling of the multiple-cut solutions according to the type of Discrete 
Fourier Transform described in [94]. 



2.2.1 Polynomials 

In this section we review the interpolation of the polynomial A(q), appearing in Eq.(p.3|), 
implemented in samurai. 

For each cut, we decompose q, namely the 4-dimensional part of g, into a specific basis 
of four massless vectors Cj [57,79,80], 



q = -po + xiei + X2e2 + xaea + ^464 , (2.6) 

such that 



ef = , ei • 63 = ei • 64 = 62 • 63 = 62 • 64 = ei • 62 = -63 • 64 , (2.7) 
and where ei and 62 are real vectors, while 63 and 64 are complex. 

The massless vectors ei and 62 can be written as a linear combination of the two external 
legs at the edges of the propagator carrying momentum q + po, say Ki and K2 , 



with 



< = ^ (i^r + ^K^) , 4 = ^ {k^2 + ^i^r ) , (2.8) 



1 - , and 7 = Ki • K2 + sgn(l, Ki • K2) J {K^ ■ K^Y - KlK^ . (2. 



2 
T 

The massless vectors 63 and 64 can be then obtained as, 

el = ^^^^ , < = ^^^^ . (2.10) 

In the case of double-cut, Ki is the momentum flowing through the corresponding 2-point 
diagram, and K2 is an arbitrary massless vector. In the case of single-cut, Ki and K2 
cannot be selected from the diagram, and are chosen as arbitrary vectors. 

After defining the basis adopted for decomposing the solutions of the multiple-cuts, 
we can list the corresponding polynomial functions, whose variables are the components of 
the loop- momentum not-constrained by the cut-conditions. 
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2.2.2 Quintuple Cut 

The residue of the quintuple-cut, Di 



Djji = 0, defined as, 

m ] 

Do--- Dn-l I 



(2.11) 



can be parametrized as [99], 



2.2.3 Quadruple Cut 

The residue of the quadruple-cut, 



^ijkeiq) = Res, 



N{q) 



%0 

= 0, defined as, 

n-l 



(2.12) 



Dn---Dr, 



(2.13) 



is parametrized as, 

+ (4:r^ + M^) [{Ks ■ e4)(« +Po) • es - {Ks ■ es){q+po) ■ e, 

=ci:r^ + cif'V' + ci:f - (4^^ + 4:f ^ /^') [(^3 • 64)^.4 - . 63)^.3 



(ei -62) , 
(2.14) 



where K3 is the third leg of the 4-point function associated to the considered quadruple-cut. 
2.2.4 Triple Cut 

The residue of the triple-cut, Di = Dj = Dk = 0, defined as. 



A ro) = Res f _ ^(gj V _ ^vkemiq)_ _\_Mq)_ 



(2.15) 



is parametrized as, 



A.,fe(g-) = 4:f +4:fM^+ 

+ 4:f ^ (9 + PO) • e4 + 4:f ' ((9 + PO) • e4)^ + 4f ^ ((9 + PO) • 64)=^ = 



+ {q + Po) • eg + 4f ^ ((9 + PO) • 63)^ + 4f ^ ((9 + f 0) • 63)^ + 



= 4:^'^ + 4:f - (4:f + 4:f ^^3) (ei • e^) + 

+ {c%'\l + c'if^t) {e, . e,f - (4f + c%'\t) (ei • e,f . 



(2.16) 
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2.2.5 Double Cut 

The residue of the double-cut, Di = Dj = 0, defined as, 

A (a) = Res [ _ ^^^) ^ V _ ^ijkjm{q)_ _ _^i_Mj.^)_ _ ^Hkio) \ ,^ 
can be interpolated by the following form, 

+ (9 + Po) • 62 + 4%^ {{q + Po) ■ 62)' + 

+ 4;5 (<? + Po) • 64 + ((g + Po) • 64)' + 

+ 4*7 ((9 +Po) • e2)((9 +Po) • 63) + ci'^siiq+Po) ■ e2)((g +Po) • 64) = 

= '^^ +C%^^'+ ('^^f -1 - ^^^^ - C%^^,) (61 • 62) + 

+ + C2axI + C%^xl - cifxiX4 - C%^XiX3^ (ci • 62)' . (2.18) 



2.2.6 Single Cut 

The residue of the single-cut, Di = 0, defined as, 

A (a) = Res I - ^^^} - V _ ^ijkjm{q)_ _ y-^ ^ijkijq) 

can be interpolated as follows, 

= 4!o + + Po) • ei) + c%{{q + po) ■ 62) + 

+ Ciliiq + po) ■ es) + Ci?4((3 + Po) ■ 64) = 

= 4!o + (cj?ia^2 + 4^)2^1 - 4^3x4 - cf^x^) (ei • 62) . (2.20) 

2.2.7 Discrete Fourier Transform 

As proposed in [94], the coefEcients of a polynomial of degree n in the variable x, say P{x), 
defined as, 

n 

P{x) = ^ Q a;^ (2.21) 

can be extracted by means of projections, according to the the Discrete Fourier Transform. 
The basic procedure is very simple: 



-8- 



1. generate the set of discrete values Pk {k = 0, 

" fc 
Pk = P{xk) = Y.^^P' e-'^'T^' , (2.22) 

by sampling P{x) at the points 

Xk = P e-'"'(^ ; (2.23) 

2. using the orthogonality relation 

N-l 

e2«A n ^-2^i\ n^^ ^ (2.24) 

ri=0 

each coefficient q finally reads, 



-^Pfce'"*(^^ (2.25) 



n + _ 

A:=0 

The extension of the DFT projection to the case of multi-variate polynomials is straight- 
forward. 

As one can notice the formula for the coefficients q, although simple, diverges when p goes 
to zero. By using the parametrization in Eq.( |2.6D , the radius p happens to be constrained 
by the on-shell cut-condition. Depending on the external invariants and internal masses, 
the dangerous value /o = might occur. In a previous work [94], we described a safer sam- 
pling, which significantly reduces the numerical instabilities arising from the vanishing of p. 
We do not repeat the same discussion here, but recall that the sampling of the multiple-cut 
solutions used for the polynomial interpolation of the triple- and double-cut residues within 
SAMURAI are chosen according to that algorithm. By using the DFT solutions as described 
in [94], we sample the numerator functions exactly as many times as the number of the 
unknown coefficients, without needing additional sampling points to improve the numeric 
precision, which would demand more computing time. 

2.3 Amplitude and Master Integrals 

The knowledge of all the coefficients appearing in the polynomials lS.ijkimi ^ijkt, ^ijk, ^ij, 
and Aj implies the following expression for the one-loop n-point amplitude, 

n-l 



. _ [d - ^) jijki) j(d+2) , (d-2)((j-4) jijH)r(rf+4)l 

i<j<k<l ^ ^ 
I ''3,0 ^ijk 2 f 

i<j ^ 

n-l 

+ j:cf>f, (2.26) 
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where, beside the scalar boxes, triangles, bubbles and tadpoles, the other master integrals 
are [84, 100] 



/ 



_ {d-2){d-^) (,+4) 

'^OT^TO- i ' ^'-''^ 

- -^r^^ , (-0) 

(£3)^ = Kff . (2.32) 

The last two master integrals, J^^^ and -f^'j-j'^ respectively a linear and a quadratic 2- 
point function, appear as a consequence of the polynomial structure of Aij(g), defined 
m Eq.(^l8|), which was chosen to have no singularity in presence of vanishing external 
invariant [89] . The vector 62 entering their definition is an element of the loop- momentum 
basis, defined in Eq.( |2.6| ), and used for the solutions of the double-cut Z)j = Dj = 0. Also, 
because of the monomial parametrization of the quintuple-cut residue, AijkemiQ), given in 
Eq.( |2.12 ), the decomposition of the amplitude in terms of MI, Eq.(2.2(;), is free of scalar 
pentagons, as already noticed in [99]. 



3. Running SAMURAI 

In this section we give some details about using samurai. All the files are available on the 
webpage: 

http : / / cern . ch/ samurai 

The archive samurai _vl . . tar . gz contains the files for the SAMURAI library, several exam- 
ples of calculations, and also the routines for the evaluation scalar integrals QCDLoop [95] 
and OneLOop [78]. 

1. Download the archive samurai_vl . . tar . gz and extract the files. They will be 
copied in a folder called /samurai. 

2. Run the Install script. It will compile all useful routines and organize them. All 
routines are written in Fortran 90 and the default compiler is gf ortran. In order 
to change compiler (or compiling options), the user should edit all the makefile 
commands. 

After running the Install script, you will find four subfolders within the /samurai 
directory: the subdirectory named /libs will contain all the libraries, namely the 
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reduction routines libsamurai . a, and three libraries for the numerical evaluation of 
the master integrals. 

Examples that reproduce all calculations described in Sec. Q can be found in separate 
subfolders in /examples. The Install script compiles all the examples, with the 
exception of the "Six Quarks" (that takes about 10 minutes to compile). The user 
can process it separately by typing make in the directory /examples/uussbb. 

3. Run each process using the corresponding command process.exe. 

The use of samurai is implemented through the following chain of calls is: 
call initsamurai (imeth, isca, verbosity , it est) 

call InitDenominators(nleg,Pi,msq,vO,mO,vl,ml, . . . , vlast ,mlast) 
call samurai (xnum , t ot , t otr , Pi , msq , nleg , rank , istop , scale2 , ok) 
call exitsamurai 

3.1 Initialization 

To initialize the samurai library, one needs to choose the arguments of the subroutine 
initsamurai 

call initsamurai (imeth, isca, verbosity , it est) 

which specify the the type of input to reduce (imeth), the routines for the numerical 
evaluation of the scalar integrals (isca), the details of the output (verbosity), and the 
test to apply to the reconstruction (it est): 

• imeth - SAMURAI can reduce integrands of one-loop amplitudes defined either as nu- 
merator of diagrams sitting on products of denominators, specified with imeth=diag; 
or as products of tree-level amplitudes sewn along cut-lines, specified with imeth=tree. 

• isca - The user can trigger the use of QCDLoop [95] by assigning isca=l, or the 
use of OneLOop [78] with isca=2. 

• verbosity - The level of information printed in the file output.dat can be chosen 
with the value of verbosity: 

verbosity=0, no output; 

verbosity=l, the coefficients are printed; 

verbosity=2, the value of the Mi's are printed as well; 

verbosity=3, the outcome of the numerical test appears. 

• it est - This option is used to select the test to monitoring the quality of the nu- 
merical reconstruction. The possibilities are: itest=0 , 1 , 2 , 3 to have respectively 
none, the global (A^ = A^)-test, the local (A^ = A^)-test, and the power-test, which 
are described in Sec.^ 

While imeth=diag supports all the options for itest, the choice imeth=tree allows 
only itest=0,2. 
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3.2 Integrand definition 

After selecting the routines for the scalar integrals and the reduction technique, the user 
should provide information about the integrand, by specifying the numerator and the 
denominators. 

The denominators of the diagram to be reduced are defined through the subroutine 
InitDenominators which generates the lists of internal momenta Pi and squared masses 
msq characterizing each propagator: 

call InitDenominators (nleg, Pi, msq, vO,mO,vl, ml, . . . , vlast ,mlast) 

The arguments of the subroutine, labeled as input/output ([i/o]) according to their role. 



are: 



nleg - [i] . The integer number of the external legs of the diagram, corresponding 
to the number of denominators. 

Pi - [o] . The array Pi(i,m) contains the nleg four-vectors present in the denomi- 



nators of the integrand, namely the vectors of Eq.(2.1) where we used the definition 
Di = {q+Pif —fnl- In the notation Pi(i,m), the first index, i=0, . . . ,nleg-l, runs 

on the set of the denominators; while the second index m=l 4, runs over the 

components of the vector, with the energy being given as 4*'^ component. 

• msq - [o] . The array msq(i), is the list of the squared masses that appear in the 
propagators. The ordering i=0 , . . . ,nleg-l is bound to the list of momenta Pi (i ,m) . 

• vO, mO - [i] . The vector vO and the mass mO are assigned to the first denominator. 

• vlast , mlast - [i] . The vector vlast and the mass mlast are assigned to the last 
denominator. 

3.3 Reduction 

Having defined the integrand denominators, characterized by Pi and msq, the actual re- 
duction of the input (xnum) is performed by the library SAMURAI, 

call samurai(xnum,tot ,totr, Pi, msq, nleg, rank, i stop, scale2, ok) 

which writes the total result of the reduction in tot. For convenience, the rational term is 
also separately written in totr. 

Here comes the detailed description of each argument: 

• xnum - [i] . The numerator of the diagram is defined in an external function, whose 
name can be decided by the user, but with fixed arguments. Hereby we adopt the 
dummy name xnum. 

The complex function xnum(icut ,q,mu2) is the integrand to be reduced. The argu- 
ments of the function xnum (i cut , q,mu2) are: 
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icut, an integer labeling the cut, where each digit corresponds to a cut-denominator 
in descending order (ex. icut= 3210 corresponds to the quadruple-cut Dq = Di = 
^2 = ^3 = ); 

q, the virtual four-momentum, q (with the energy given as 4*^^ component); 
and mu2 the extra-dimensional mass-scale, ^j?. 

When inieth=diag, xnum is expected to have the form of a numerator, hence being 
polynomial in q and ji^. In this case xnum is a unique function to be processed at 
every level of the top-down reduction by cycling on icut, but does not depend on 
the considered cut. 

When imeth=tree, xnum is expected to be formed by the product of tree-amplitudes, 
therefore the presence of propagators it is also allowed. In this case, xnum is not 
unique, but should change according to the considered cut. Therefore, the value of 
icut yields a selective access to the proper integrand within the same function. 

• tot - [o] . The complex variable tot contains the final result for the integrated 
amplitude of numerator xnum. The finite part, that also includes the rational term, 
will be stored in tot(O), while tot(-l) and tot (-2) contain the single and double 
poles, respectively. 

• totr - [o] . For the purpose of comparisons and debugging, we also provide the ra- 
tional part totr alone. This complex number is the sum of all contributions coming 
from integrals in shifted dimensions, namely all contributions that contain a depen- 
dence from /x^ in the reconstructed integrand. 

• nleg - [i] . Already defined. 

• Pi - [i] . Already defined. 

• msq - [i] . Already defined. 

• rank - [i] . This integer value is the maximum rank of the numerator. This in- 
formation is extremely valuable in order to optimize the reduction and improve the 
stability of the results. Using this information, we can simplify the reconstruction of 

the numerator by eliminating contributions that do not appear in the reduction. If 
the information about the rank of the integrand is not available, rank should be set 
equal to nleg. 

• istop - [i] . This flag stops the reduction at the level requested by the user, istop 
is an integer, whose range of values is from 1 to 5. 

istop=5,4,3,2,l will interrupt the calculation after determining pentagon, box, 
triangle, bubble, and tadpole coefficients respectively. This procedure can be partic- 
ularly useful to improve the precision of calculations when one knows a priori that a 
particular set of integrals does not contribute to the considered process. 
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• scale2 - [i] . This is the scale (squared) that is used in the evaluation of scalar 
integrals. 

• ok - [o] . This logical variable carries information about the goodness of the re- 
construction. The default values is ok=true, and it is set to ok=false when the 
reconstrucion test fails. 



As stated in Section 4.7.1, the generic one-loop integrand can be polynomial in e up 
to the second-order. Each coefficient of the e-decomposition can be assigned to a specific 
function, i.e. xnumO, xnuml, xnuni2, which can be independently processed. 

3.4 Reconstruction Tests 

There are three different ways of monitoring the quality of the coefficients reconstructed 
by SAMURAI. 

3.4.1 Global {N = iV)-test 

The first option (itest=l) is the so-called 'W = A^" test on the reconstructed expression for 
the numerator functions, which was already discussed in [80,89]. It is based on the equality 



given by Eq.(2.4), between the original numerator in the l.h.s. and the reconstructed one 
in the r.h.s., evaluated at an arbitrary value of q. 

A possible drawback of this precision test lies in the fact that the coefficients of tadpoles 



and bubbles in Eq.(2.4) multiply a large set of denominators: for a six-point function, each 
tadpole coefficient multiplies five denominators, namely a term proportional to masses or 
momenta, q, raised to ten powers, that can be huge in some cases or very small in other 
situations. This might have the effect of hiding the contribution of some coefficients or, as 
happens more frequently, might yield to overestimating the error in the reconstruction. 

3.4.2 Local {N = A)-test 

A second check is a "local N = N" test (itest=2), regarding the reconstruction of each 
polynomial A(g), respectively defined in Eqs.(2.11, p. 13 , 2.15, p. 17 , 2.19| ). In this case 



the value of q used for the numerical check is chosen among additional solutions of the 
considered multiple-cut, which have not participated to the determination of A(g) itself. 
This option is suitable for a unitarity-based calculation (imeth=tree). 

3.4.3 Power-test 

A third option (itest=3) for testing the precision of the reconstruction is the "power test". 
We can observe that the maximum powers in q in the r.h.s and l.h.s of Eq.( |2.4D are differ- 
ent: the reconstructed side can contain terms with high powers of q that are not present in 
the original numerator. Therefore it is clear that the overall coefficients in front of these 
terms should vanish. 

The reconstructed expressions in general are not simple, since they involve pieces coming 
from the polynomial spurious terms multiplied by the denominators. However, for each 
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choice of the rank and number of denominators, there is at least one simple set of coeffi- 
cients that sum to zero exactly. Moreover, this set is the lowest one in the reconstruction 
and therefore it carries information about any loss of precision at previous steps in the 
reduction. 

If the difference between the rank and the number of denominator is equal to three 
(nleg-irank =3), the sum of all the coefficients of three-point scalar integrals should be 
zero, namely: 

E^3.(0) = (3.1) 

i 

where the sum is over all possible triple cuts. 

Analogously, if the difference between the rank and number of denominator is equal to two 
(nleg-irank=2), the sum of the coefficients of two-point scalar integrals should be zero, 
namely: 

E^2,^(0) = (3.2) 

i 

where the sum involves all double cuts. 

Finally, if the difference between the rank and number of denominator is equal to one, 
(nleg-irank=l), the sum of the coefficients of the tadpole scalar integrals should be zero, 
namely: 

^ci,i(0) = (3.3) 

i 

where the sum involves all single cuts. 

The situation is slightly more complicated for maximum rank when difference between the 
rank and number of denominator is equal to zero. If (nleg-irank=0), we should consider 
all the one-point spurious coefficients ci^j(l) to ci^j(4), each multiplied by the corresponding 
vector ei^j, . . . , e^^i of the basis defined in Section |2.2.1 . Summing over all possible single 
cuts, labeled by i, we get the condition 

EE^mH<. = (3.4) 

i n=l 

As a final remark, we observe that the outcome of the "power test" does not depend in 
any way from the choice of the integrated momentum q, unlike the previous two methods. 

The threshold values for the reconstruction checks can be set in the file Itest .dat, to be 
located in the directory where the call to initsamurai is made. The phase-space points 
failing the tests (ok=false) are stored in the file bad. points, in the same directory. In 
principle they could be re-processed enhancing the numerical precision by compiling the 
SAMURAI library in quadruple-precision. 

3.5 Comments on Precision 

The precision of the results obtained using a reduction algorithm at the integrand-level 
depends on many variables. 
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When the numerator is a real function of the external momenta and masses there 
is a simple way to establish the quality of the reduction: real functions give rise to real 
coefficients of MI. In this case, the error on each coefficient can be estimated by the size of 
the imaginary part, that should vanish. 

More generally, the quality of the reconstruction can be quantified by the ratio of 
the difference between the exact calculation (analytical or multi-precision) and the recon- 
structed one, and the former, evaluated over a large set of unweighted points. This proce- 
dure gives a good indication, but it is not always safe, because the error on the prediction 
in a calculation based on the importance sampling could suffer from the accumulation of 
bad points in the neighborhood of higher weights. 

We identify three kinds of possible instabilities, which could be all controlled by adopt- 
ing quadrupole or multiple precision routines. 

The first kind of instabilities is related to the well known problem of the vanishing of 
the Gram determinants, inducing an enhancement of the coefficients of the MI carrying 
such pathological kinematic factor. They can be monitored by the tests implemented in 
SAMURAI, and the dangerous cases could be dealt with by introducing branches to dedicated 
reduction routines, hence without making use of the multiple precision. 
The second kind corresponds to big cancellations among the contributions from different di- 
agrams in the same calculation. On-shell methods, which work with purely gauge invariant 
objects, seems to represent the best option to avoid such problem. 

The third type of instability can occur when the values of internal masses are sensibly 
larger then the phase-space invariants. In this case, both the cut-constructible part and 
the rational term are large but their sum remains relatively small. This in principle could 
be cured with a change of the integral basis where the cancellations are built-in. 

Our tool does not switch automatically between double and quadruple precision. The 
running in the latter case is time-consuming, therefore, along the lines of the above consider- 
ations, we are investigating a more systematic treatment of the problematic configurations, 
which goes beyond the scope of this version of the code, and will be the subject of a future 
publication. 

4. Examples of Applications 

In this section we present examples of calculations of one-loop amplitudes performed with 
SAMURAI. These examples axe selected with the idea of covering different situations and 
problems that can be treated within our code. Our intention is to show the fiexibility of 
this framework and present examples of applications performed in various regularization 
schemes widely used for the calculation of one loop virtual corrections. 

SAMURAI can process two different kinds of input, according to the strategy adopted 
for the generation of the integrand. In the Feynman diagrams approach one should provide 
a set of numerator functions, each accompanied by a corresponding list of denominators. 
On the other hand, in the generalized unitarity approach the input will be in the form of 
products of tree-level amplitudes. In the following we describe some calculations performed 
within both frameworks. 
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In several cases, we use Rambo [101] for generating phase-space points. 



4.1 Four-photon Amplitudes 

This example is useful to verify the proper reconstruction of the rational term. The leading 
term of the process 77 — )• 77 in QED proceeds via fermion-loop [102, 103]. We treat both 
cases of massless and massive fermion. The four-photon amplitudes get contributions 
from the 6 Feynman diagrams representing the possible permutations of the 4 photons 
attached to the fermion loop. Indeed, only 3 permutations are independent and need to 
be evaluated, because loops related by flipping the fermion line give the same answer. Let 
us consider the diagram with the photons labeled in clockwise order 1234, carrying the 
following denominators, 

{Lj - m^) m - m^) (Li - rr?) {Lj - m^) (4.1) 

and numerator, 

Niq) = -Tr [(^^ + m) /a ih + ^) h ih + ^) h (^4 + ^) ^1] (4-2) 

where 

Li = q , L2 = q+P2 , L3 = q + P23 , L4 = q + p234 , (4.3) 

with pjjfc = Pi+Pj+Pk- The other two independent contributions are obtained by permuting 
momenta and polarizations: (234) — t- (243), (324). So, in this example, the inputs to run 
SAMURAI are simply: 'Pi = {0,P2,P2 + P3tP2 + P3 + Pi), msq = (m^, m-^, m^, m^), irank = 
4, istop = 1. Once the loop momentum is decomposed as in Eq.( p.2[ ), we end up with an 
expression suitable for the numerical evaluation: 

Niq, /i^) = -(m^ -fi^m^ + /) Tr[/2 k h h\ 

- (vr? - iJ?) (TY\i2 h h U h h] + Tr[/2 h h h h h\ 
+ Tr[/2 h As h A ii] + Tr[/2 f^2 h h h Ai] 

+ Tr[/2 1^2 h ii A /i] + Tr[/2 J^2 h h A i\\ ) 

- Tr[A/2A/3A/4A/i] , (4.4) 

where Lj is the 4-dimensional part of Lj (= Li + fi). Note that now the whole expression 
can be evaluated numerically in terms of the four dimensional complex variable q and the 
real variable /i^. Using SAMURAi, it is easy to see that: the term proportional to /^^ 
in Eq.( |4.4D gives rise to null integrals and does not contribute, the terms proportional to 
^2gii g^j^g j^Q^ individually zero but they cancel when summing over all contributions; and 
finally that the fi'^-teim gives the correct rational term. 

4.2 Six-photon Amplitudes 

The six-photon amplitudes [83, 89, 104-108] are also a good test for the reconstruction of 
the rational term, that, after summing over all diagrams, has to vanish [109]. 
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The construction of the amphtudcs follow closely the one that we used for the four photons. 
Out of the 120 contributing diagrams, all containing up to rank-6 tensor integrals, only 60 
need to be computed. We can construct all of them as permutation of just one diagram. In 
the massless case, we consider the diagram with the photons in the clockwise order 123456, 
whose corresponding numerator reads, 

N{q, fi^) = -Tr [f., ^ U h h h U h h h U A] • (4-5) 

where 

Li = q , L2 = q + P2 , L3 = q + P23 , Li = q + p2u , 

L5= q+ P2345 , Lq = q+ P23456 • (4.6) 

This example turns out to be challenging for the reduction algorithm, because each diagram 
separately admits a non-trivial reduction with non- vanishing coefficients for all the MI and 
rational terms but, after summing together the partial results of all diagrams, there are 
strong cancellations. In the final answer all contributions coming from 2-point functions 
cancel out. Moreover, also the rational terms vanish. 

Indeed, the final expression contains only cut-constructible terms and no rational part and 
the knowledge of the coefficients of boxes and triangles alone is sufficient to obtain the 
correct answer for the total amplitude. 

After the dimensional decomposition of the loop momentum q, it is easy to see that 
all the terms containing one, two or three powers of /i^ give rise to vanishing integrals and 
do not contribute. As a consequence, the only term needed in the numerical evaluation is 
the four dimensional one: 

N{q, /i^) = N{q) = -Tr [f., ^2 h h h h U h h h U ■ (4-7) 

For a numerical check we consider the value of the two amplitudes ^(— ,— ,+, +,+,+) and 
^(+, -,-,+,+,-) [107,108] 

4t -, +, +, +, +) = 11075.04009210435 , (4.8) 

a-' 

4t ^(+, -, -, +, +, -) = 7814.762085902767 , (4.9) 
ay 

evaluated at the phase-space point [104], 

P3 = (33.5,15.9,25.0) (4.10) 

P4 = (-12.5,15.3,0.3) (4.11) 

P5 = (-10.0,-18.0,-3.3) (4.12) 

P6 = (-11.0,-13.2,-22.0) (4.13) 

with px and p2 directed along the positive and negative z-axis respectively. 

By running SAMURAI with istop = 2, namely keeping the contributions of the bubbles, 
the results are: 

^ A{-, -, +, +, +, +) = 11075.040 174990 , (4.14) 
^ ^(+, -, -, +, +, -) = 7814.762 3429908 . (4.15) 



ar 
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Figure 1: Results for the 8-photon amplitude with helicity h + + + ++. The continuous line 

represents the analytic results of [83] . The results of SAMURAI are produced in double-precision and 
with istop=4. 

By exploiting the knowledge that contributions from bubbles and rational terms will vanish, 
and therefore removing these terms from the reduction, we verify an improvement on the 
final result. Infact, by setting istop = 3 and isolating only the cut-constructible terms (by 
subtracting totr diagram by diagram), the results of samurai turn out to be in better 
agreement: 

-4 -.+,+,+■+) = 11075.04009210 2 . (4.16) 

-4 -■-■+.+,-') = 7814.76208590 84 . (4.17) 

As expected, the strong cancellations between the 60 diagrams spoil the precision of 
the full results even if the number of good digits for this specific phase-space point can still 
be considered sufficient for phenomenological studies. 

4.3 Eight-photon Amphtudes 

The eight-photon amplitudes [83,106,110] are an example of the functionality of SAMURAI 
for many-particle scattering. 

The numerator function is written along the same lines as in the previous two sections. 
In this case, the number of diagrams is 5040. We evaluate the amplitudes for two helicity 
choices. 

By using the same sampling set as in [106], we show in Fig.|l] how the numerical result 

produced with samurai in the MHV case, h + + + ++, are tight to the analytic 

behavior [83]. The NNMHV case, h + ++, shown in Figj2|, is a new result that 

confirms the structure of the amplitude discussed in [110], where only boxes do contribute. 
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Figure 2: Results for the 8-photon amplitude with helicity h +++, produced with samurai 

in double-precision and with istop=4. 21 points are also given with samurai in quadruple-precision 
and with istop=2. 



4.4 Drell-Yan 



Figure 3: Triangle diagram for Drell-Yan. 



The one-loop correction to uu — )• e~^e [111, 112] is an easy example of a numerator 
with e-dependent terms. The numerator of the diagram in Fig. ^ can be cast in the form 

N{q, = Cf gl n(Pe- ) ^ v{Pe+) HPu) [ 2 (2 - d) | + [ (d - 2) 
+4:ipu-q-Pu-q-Pu- Pu) ] 7^ ] u{pu) 

with denominators 

\2 



(Q+Pu) {q + Pu+Pe-+P, 



The value d = 4 in the expression above corresponds to the result in the Dimensional 
Reduction (DR) scheme, while the choice d = 4 — 2e yields an e-dependent term, according 
to the Conventional Dimensional Regularization (CDR) scheme, samurai can be used 
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Figure 4: Parent diagram for the amplitudes A^-i. 



to reduce both the e and the coefficient of the e term individuahy, namely Nq and A^i 
of Eq.(p.3|). It is easy to see that the inclusion of the latter has the well known effect of 
subtracting a contribution Cp gl times the tree- level amplitude from the finite part of the 
DR-result. 

4.5 Leading-color Amplitude for V + Ijet 

The leading color amplitude for the virtual NLO correction to ^ + Ijet production at the 
hadron collider is a good exercise to show the reduction in a case where the contribution 
of all diagrams is cast in a single numerator function. 

Once the color factors have been stripped, this amplitude can be calculated at the 
Feynman diagram level taking the sum of the parent diagram in Fig. ^ and its pinched 
diagrams, i.e. four triangles and two bubbles. The presence of the 7^ in the weak vertex 
imposes a choice on its treatment in dimensional regularization. Adopting the Dimensional 
Reduction (DR) scheme and assuming an anticommuting 7^ one can get the right result 
adding a well known finite-renormalization contribution, amounting to {—Nc/2) times the 
tree-level amplitude. 

With the proper routing of the loop momentum in the diagrams, it is possible to collect 
all the diagrams over the four denominators of the parent box: the numerator of triangles 
is multiplied by the single missing denominator, while the bubbles by two denominators. 
In this, we should process only one numerator function. This way of collecting the dia- 
grams does not spoil the precision of the result. Using this construction, we found perfect 
agreement with the expression for ^5;! given in the Eqs.(D.l-D.5) of [113]. 

4.6 Five- and Six-gluon amplitudes 

We choose two simple examples, namely the amplitudes contributing to the rational part of 
the all-plus helicity 5-gluon and 6-gluon scattering [60,114-116], to show how a unitarity- 
based calculation can be implemented within samurai (option iineth=tree). 

The diagrams involved correspond to one-loop amplitudes with external gluons coupled 
to a massive-scalar loop, whose integrand can be built by means of the tree- level amplitudes 
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given in [116,117], namely 



4-^(1,; 2+; 3,) = !^^ , (4.18) 

^tree.i .3+ 3+ 4+. 5 /.M2| 1 (2 + 3) |4] 

^ ^ ^' ' ' ''^^-(23)(34)(pf,-M2)(pl,-;.2) ' (4-20) 

where r2 is the reference vector of the gluon-2, and pij = ki + kj. For instance, the integrand 
of the quintuple-cut shown in Fig.|5| can be written as, 

N{q, fi^)=Ai{Lr, 1+, 2+; -L2) x ^3(^2; 3+; -L3) x ^slLs; 4+; -L4) 

x^3(^4; 5+; -L5) X ^3(^5; 6+; -Li) (4.21) 

where 

Li = q , L2 = q + pi2 , L3 = q + P123 , = q + P1234 , L5 = q - pQ . (4.22) 

In the case of the 5-gluon amplitudes, we give the complete set of integrands, for quintuple-, 
quadruple-, triple- and double-cuts (istop=2), although only boxes appear in the result. 
In this case, we see explicitly that triangles and bubbles do not contribute. For the same 




Figure 5: Quintuple-cut of the 6-gluon amplitude. 



reason, in the 6-gluons case, we only give the integrands for the quintuple- and quadruple- 
cuts (istop=4). 

The results of these calculations, due to the external helicity choice, are purely rational 
in the d = 4 limit and agrees with the results of [117]. 

4.7 Six-Quarks Scattering 

When the number of diagrams contributing to the scattering amplitude is small, the input 
file that includes the numerators and the list of momenta to be processed by the reduc- 
tion is fairly simple, and the calculations are feasible with a minimal amount of automa- 
tion [118,119]. However, even in simple careful automation reduces the probability 
of introducing bugs or human mistakes in the code. 
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An automatized generation of the input files becomes a necessity as tlie complexity of 
the process increases. As a final example (with the diagrammatic approach), we tackle a 
more involved calculation, namely the one-loop QCD corrections to the 6-quark scattering 
qiqi —7- q2Q2Q3Q3- The number of Feynman diagrams contributing to this process requires 
a fully automated approach. 

The amplitude for qiqi — )• q2q2Q3Q3 involves 258 Feynman diagrams (8 hexagons, 24 pen- 
tagons, 42 boxes, 70 triangles, and 114 bubbles). Each diagram, or convenient combinations 
of them, should be processed by the reduction algorithm separately. The numerators and 
the lists of denominators required by the reduction have been generated and automatically 
written in a Fortran90 code, ready to be processed by samurai. 

We use this example also as a first benchmark on the functionality of our framework. 
During the generation of the code, all Feynman diagrams contributing to the process are 
automatically written and organized in Fortran90 files fully compatible with the reduction 
library, ready to be run. 

In order to check our algebraic manipulations, we compute both NQ{q) and Ni{q) of 
Eq.(p.3|), namely also the part of the numerator proportional to e, although in an actual 
calculation this can be avoided by choosing the regularization scheme conveniently. 

There are eight different helicity configurations that contribute to this process. Our 
numerical results have been compared with those obtained for the same process with 
golein-2 . and golem95 [96] and we found perfect agreement. 

On a Intel(R) Xeon(R) CPU X5482 3.20GHz machine, the generation of the code for 
the full process takes less than 10 minutes, and the result for each color summed helicity 
amplitude is produced in 55 ms per phase-space point. However, by avoiding the reduction 
of Ni{q) with a proper scheme choice, the computing time goes down to 36 ms/ps-point. 

4.7.1 Numerator 

When working with Feynman diagrams, we prepare the numerator function A^(^) by pro- 
cessing the output of a diagram generator symbolically with a computer algebra program; 
the actual computer program is written by an optimizing code generator (see also Fig. |6|). 
This modular approach is very generic and, to a large extend, can be based on existing tools; 
in particular we have an automated setup using QGraf [120], Form [121] and haggles [122]. 
Furthermore, the matrix element generator golem-2.0 [96] has been extended to provide 
an interface which simplifies the use of the components mentioned above. We want to 
stress that the described setup is very modular and that any component in the workflow 
can be exchanged by alternative solutions. 

As discussed in Section ^, the most general numerator of one-loop amplitudes, J\f{q, e), 
can be written as, 

M{q, e) = No{q) + eNi{q) + e^N2{q). (4.23) 

The functions A^'o, Ni and N2 are functions of q'^ and fi^, therefore in our discussion, except 
when it is necessarily required a distinction between them, we will simply talk about N, 
giving as understood that the same logic would apply to each of the three contributions N^. 
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Figure 6: Schematic overview of the generation of the numerator. The boxes correspond to the 
steps required for the generation of the numerator for a given process using Feynman diagrams. 
The actual implementations we have used for each step are given in rounded rectangles (see text). 
The dashed arrows indicate control flow, plain arrows indicate data flow. 

We work with the hclicity projections of the amplitude which are decomposed into 
subamplitudes formed by the sum of all diagrams sharing the same set of denominators. 
The color information is hidden from the reduction by defining the numerators of the 
subamplitudes, J\f^'^\q,e), as the contraction of the numerators of the one-loop diagrams 
with the tree-level amplitude. If we call N^^^^^ the numerator stemming from the sum 
of all diagrams which have (exactly) the denominator Di^ Di^ ■ ■ ■ Di^ , the corresponding 
subamplitude would be 

^^'\^^^)=^Ln■N{^,^,...^^}. (4.24) 

In our implementation this product is done numerically and does therefore not add to the 
complexity of the expressions. In cases where the tree- level matrix element vanishes, one 
can always find an appropriate set of color projectors v\Vi into one-dimensional subspaces 
such that 

A■A^ = Y.^VlAn)'^ ■ {VlAn). (4.25) 

I 

where Vi correspond to Wigner-Eckhard symbols. In the cases with no external color, the 
only projection is Vq = 1. The objects Vi ■ -/Vjj^jj j^j hence are the objects that undergo 
the reduction. 

Optionally, one can also group larger sets of diagrams into subamplitudes by also 
considering diagrams which contain a subset of the maximal set of denominators. The 
numerator of the corresponding subamplitude in the latter sense would be 

= [N{hi2...ir.} + DinN{i,i,...i„_,y + A„_iA^{n...i„_,i„} + ■ ■ ■ 

+Ai • • • A„_.-DniV{,„_,} + Ai A. ■ ■ ■ A„_iiV{.„}] • 4orn- (4-26) 
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4.7.2 Algebraic Simplification of the Lorentz Structure 

In order to unravel the dependence of M{q, e) on q, fj,"^ and e we use dimension sphtting 
based on the 't Hooft-Veltman scheme. We define the subspaces of the regulated Minkowski 
space such that 

r'' = 9'"' + r'', 9^1 = d, g^^ = A, ~9ll = -2e, 5^^%. = (4.27) 

and with the corresponding projections of the Dirac matrices 7^ = gu^'^ and 7^ = Quj'^ 
the Dirac algebra is uniquely defined by 

{7^r} = 25'^^ {^M,^5} = 0, [7^75] = 0, (4.28) 

Working with this scheme one can show [123] that after separating the four from the 
{d — 4) dimensional projection of each Dirac matrix one can factorize a mixed spinor line 
into 



(pa|7''' •••7'"=7'"=+' = (pa|7^' •••7'"=|Pa') • tr{7'''=+i •••7'''}/tr{l} . (4.29) 

In this notation the definition of the helicities is such that \p±) = |(1 ± "f5)u{p) and 
{p±\ = 'i*(p)^(l =t 75)) where p and p' are lightlike vectors. The extension to massive 
vectors is straightforward by projecting each massive vector onto a sum of two lightlike 
vectors. The trace in Eq. ( 4.29| ) evaluates to a product of metric tensors g^'^J using the 



usual rules for spinor traces. Since in the 't Hooft-Veltman scheme at one-loop the only 
d-dimensional vector is the integration momentum these metric tensors lead to factors 
of fi^ and e. The Lorentz indices inside the remaining, four-dimensional spinor lines are 
eliminated using Chisholm identities, of which we apply also a variant specific to spinor 
chains, 

{px\T^^T'\p',,) .J, = 2 (T'\p',,){p,\r - XX''V\px){p'y\'V') (4.30) 



where F and V are strings of four-dimensional Dirac matrices and ^ denotes the string in 
reversed order. 

After these steps, the numerator is suitable for efficient numerical evaluation since it is 
expressed entirely in terms of constants, dot products and spinor products of the form 
{Px\p'x) and {px\4\p'-x)- 

4.7.3 Result of ^(1", 2+, 3", 4+, 5", 6+) 



The LO contribution and the NLO virtual corrections to the squared amplitude (ultra- 
violet renormalised) are defined as, 

OLO = -^LO-^LO (4.31) 

A . r Ois i^T^Y (0,-2 fl-l \ N 
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The result of goleni95 for the hehcity configuration {q^ ,Qi,Q2 ^it ^% ^lt)^ ps-point 
given in Eq.( [4.13| ), is 

OLo = 0.9686295685264447 x 10"^ , (4.33) 

a_2 = -8.000000000048633 , (4.34) 

a_i = 46.40675046335535 , (4.35) 

ao = -233.8908276457752 ; (4.36) 

and the one computed by samurai is 

flLo = 0.9686295685264458 x 10"^ , (4.37) 

a_2 = -7.999999999999935 , (4.38) 

a_i = 46.40675045992446 , (4.39) 

ao = -233.8908276128404 , (4.40) 

showing a nice agreement (the color-average factor, 1/9, and the helicity- average factor, 
1/4, are already included). 

The double- and single-pole of the virtual contribution are consistent with the infrared 
poles amounting to [28], 

a_2 = 8.000000000000000 , (4.41) 
a_i = -46.40675046319159 . (4.42) 

4.7.4 Precision of Integrated Results 

We have used the matrix element of the qiqi — )• q2Q2Q3Q3 amplitude for recalculating the 
qiqi — )• ^2^292^2 amplitude [96] by anti-symmetrizing over the final state. We have in- 
tegrated the virtual matrix element with MadEvent [2, 3] and compared the poles of the 
virtual amplitude to those of the integrated dipoles using MadDipole [34,36]. Figure ^ 
shows the remainder of the pole contributions which should sum up to zero. The results 
represent a realistic Monte Carlo integration and indicate that the precision is well under 
control. 



5. Conclusions 

In this work we have presented samurai, a tool for the automated numerical evaluation 
of one-loop corrections to any scattering amplitudes within the dimensional-regularization 
scheme. Its implementation is based on the decomposition of the integrand according to the 
OPP-approach, extended to the framework of the generalized d-dimensional unitarity-cuts 
technique, and on the use of the Discrete Fourier Transform as polynomial interpolation 
technique. We have shown how samurai can process integrands written either as nu- 
merator of Feynman integrals, like in diagrammatic methods, or as product of tree-level 
amplitudes, according to unitarity-based methods. In both cases, the advantage of working 
within a d-dimensional unitarity framework is that the result of SAMURAi is complete and 
does not require any additional information for the reconstruction of the rational terms. 
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Figure 7: Estimate for the precision obtained from the difference between the single (resp. double) 

poles of the virtual amplitude and those of the integrated dipoles for qiqi — > (72 92 92 92- The results 
have been obtained by integrating 10^ phase-space points at y/s = 14TcV, where we have used 
cuts on pt > 30GeV and the rapidity r] < 2.5 as well as a separation cut of AR > 0.8 between 
the final state particles. We used the CTEQ6m [124] PDF set with two-loop running for with a 
renormalisation scale of n = \/^iPT{iY ■ 

We discussed its application to a series of examples such as the 4-, 6-, and 8-photon 
scattering amplitudes in QED, the QCD virtual corrections to Drell-Yan, the leading 
color amplitude for V + Ijet production, the six-quark amplitudes, and contributions from 
massive-scalar loop to the all-plus helicity 5- and the 6-gluon amplitudes. For the six- 
quark scattering gigi — > 52^293^3, we also considered a fully automated reduction, from the 
integrand generation to the final result. 

Given the versatility of the code, SAMURAI may constitute a useful module for the 
systematic evaluation of the virtual corrections, oriented towards the automation of next- 
to-leading order calculations relevant for the LHC phenomenology. 

The reduction library libsamurai and the examples are publicly available at the web- 
page: 

http : / / cern . ch/ samurai 
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